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Abstract 

In this paper we present qualitative and quantitative comparison 
of various analytical and numerical approximation methods for cal- 
culating a position of the early exercise boundary of the American 
put option paying zero dividends. First we analyze their asymptotic 
behavior close to expiration. In the second part of the paper, we intro- 
duce a new numerical scheme for computing the entire early exercise 
boundary. The local iterative numerical scheme is based on a solu- 
tion to a nonlinear integral equation. We compare numerical results 
obtained by the new method to those of the projected successive over 
relaxation method and the analytical approximation formula recently 
derived by Zhu. 
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1 Introduction 

The analysis of the early exercise boundary and the optimal stopping time 
for American put options on assets paying zero dividends has attracted a 
lot of attention from both theoretical as well as practical point of view. 
An American put option is a financial contract between the writer and the 
holder of the option. It gives the holder the right, but not the obligation, 
to sell the underlying asset at the prescribed strike price any time before 
expiration. Under the standard assumptions made on the underlying stock 
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process and completeness of the financial market (c.f. [T81I23] ) the American 
put option can be priced using the Black-Scholes equation (c.f. [1]) on a time 
dependent domain of the underlying asset price. More precisely, the early 
exercise boundary problem for the American put option can be formulated as 
follows: find a solution V = V(S, t) and the early exercise boundary position 
Sf = Sf(t) satisfying 

~dt +rS ds + T S ds 2 ~~ rV = ' 0<t<T ' S f (t)<S<oo, 
V(+oo,t) = 0, V{S f {t),t) = E-S f (t), ^(S f (t),t) = -1, (1) 



V(S, T) = (E-sy 



The solution V(S,t) is defined on a time-dependent domain S G (Sf(t), oo), 
where t G (0, T) (cf. Kwok [23]). Here S > stands for the underlying stock 
price, E > is the exercise (strike) price, r > is the risk-free rate, a > 
is the volatility of the underlying stock process and T denotes the time of 
maturity. In what follows, we denote by r = T — t the time to maturity. 
The function [0, T] 3 t t— > Sf(t) G M represents the early exercise boundary 
position. The above mathematical formulation of the problem of pricing the 
American put option by means of a solution to the free boundary problem is 
a basis for development of various integral equations for describing the early 
exercise boundary position Sf(t). The analytical approximation formulae 
are often based on approximation of a solution to such an integral equation. 
Notice that there are also other numerical methods for approaching the free 
boundary problem (j2J) like e.g. front-fixing and transformation methods. We 
refer the reader to papers by Kwok and Wu [21], Sevcovic [221 ED|, Ankudi- 
nova and Ehrhardt [2] and references therein. Following Kwok [23], a solution 
V = V(S,t) to the problem of pricing the American put option fulfills the 
following variational inequality: 

(W F)V a 2 8 2 V 

-W +rS gS + Y i 'd¥~ rV -°- v (S,t)>V(S,T), 

for all < t < T, < S < oo, (2) 
V(0, t) = E, V{+oo, t) = 0, for < t < T , 
V(S, T) = {E- S) + , for < S < oo. 

The formulation of the problem of pricing American put option as a varia- 
tional inequality is often used when we need to compute not only the free 
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boundary position Sf(t) but also the entire solution V(S, t). The above varia- 
tional inequality can be effectively solved by means of the so-called projected 
successive over relaxation (PSOR) method by Elliot and Ockendon [T3j. 

In the last decades, many different, but equivalent, integral equations for 
pricing the American put option have been derived by Barone-Adesi and 
Whaley [3], Bunch and Johnson [5], Carr, Jarrow and Mynemi [5|, MacMil- 
lan [2B] and others. The asymptotic analysis often leads to an approximate 
expression of the free boundary close to expiry. Since the closed form ana- 
lytical formula for the early exercise boundary position is not known, many 
authors (see e.g. Geske, Johnson and Roll [El EE], Johnson [19], Karatzas 
|21j . Evans, Kuske and Keller [221 EH], Mynemi [28] and recent papers by 
Alobaidi et al. [UEZ], Stamicar et al. [31] . the survey paper by Chadam [7] 
and other references therein) investigated various approximation models and 
derived different approximate expressions for valuing American call and put 
options. We also refer to the books by Kwok [23] and Wilmott et al. [32] 
for a survey of classical theoretical and computational results in the field of 
pricing the American put option. 

In this paper, we focus on comparison of the valuation formulae due to 
Evans, Kuske and Keller [221 [15], Stamicar, Sevcovic and Chadam [31J and 
the recent analytic approximation formula by Zhu [33] (see also [3^1135] ). Our 
main goal is to present qualitative and quantitative comparison of the above 
mentioned analytical and numerical approximation methods for calculating 
the early exercise boundary position. In the first part of the paper, we analyze 
and compare asymptotic behavior of the early exercise boundary close to 
expiry for analytical approximations developed by Stamicar, Sevcovic and 
Chadam [31], Evans, Kuske, Keller [221 US] and Zhu [331 EH 135] • We show 
that the approximation formulae due to Evans, Kuske and Keller [2"2"1 [T5] 
and Stamicar, Sevcovic and Chadam [3T] have the same asymptotic behavior 
of Sf(t) as t — > T. We also show that the analytic approximation formula 
due to Zhu has an asymptotic behavior differing from the previous ones by 
a logarithmic factor. In the second part we propose a new numerical scheme 
for computation of the entire function Sf(T),t G [0,T], based on a solution 
to the nonlinear integral equation from [21]. We compare numerical results 
obtained by the new numerical method to those of the projected successive 
over relaxation method by Elliot and Ockendon [H] for solving the variational 
inequality (j2J) and the analytical approximation formula recently developed 
by Zhu et al. in [331 EU [35] . 
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2 Analytical approximate valuation formulae 



In this section we present a survey of analytical, implicit integral and nu- 
merical approximation schemes for computing the early exercise boundary 
for the American put option. First we focus on the recent result due to Zhu 
who in [33] derived a closed analytic approximation formula for the early 
exercise boundary position Sf(t) = g(T — t). We also derive the asymptotic 
behavior of Sf(t) for t — > T. Next we concentrate on implicit representation 
formulae for g(r) expressed in the form of a single nonlinear integral equation 
for the function g. We recall implicit integral equation derived by Stamicar, 
Sevcovic and Chadam in [3T]. We again derive the asymptotic behavior of 
the early exercise boundary position as t — > T. In the last subsection we 
present another approximations derived by Evans, Kuske and Keller [22l [15] . 



2.1 Analytical approximation valuation formula by Zhu 



In this section we recall a recent interesting result due to Zhu. In [33J Zhu de- 
rived a new analytical approximation formula of the early exercise boundary 
by application of the Laplace and inverse Laplace integral transforms to a 
dimensionless form of the governing parabolic PDE and successfully obtained 
a closed analytic approximation formula for the early exercise boundary posi- 
tion as a sum of a perpetual option and integral that valuates early exercise 
boundary position. The resulting formula for the early exercise boundary 
Sf(t) = g(T — t) reads as follows: 



g Zhu (r) 



~iE 2E f°° Ce~ T 



7 



TT 



-(a 2 +C 2 ) 



a 2 + C 2 



-e-^ S in(/ 2 *(C)R, (3) 



where 7 



1+7 
2 ' 



1-7 
2 ' 



and 



mo 
/2(c) 



b 2 + c 2 

1 

6 2 + C 2 



b In ( -Va 2 + C 



Cm -Va 2 + C 2 

7 



( arctan(C/a) 



6arctan(C/a) 



(4) 



Notice that the first summand in ([3]) represents the constant value of a per- 
petual put option i.e. the limit lim T ^ 00 g(r) = jE/(l + 7). 
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Early exercise boundary asymptotic close to expiry 

Next we examine the asymptotic behavior of the function g Zhu {r) for r — > 0. 
Notice that we have g(0) = Sf(T) = E (c.f. Kwok [23]). We shall prove that 

E-g Zhu (r) 1 „ 
lim V / = -=Ea. 

r^o+ V^-lnr) V27T 

Indeed, if we introduce the change of variables: s = r £ ^-(a 2 + £ 2 ) we obtain 

E - g Zhu (r) 2E [°° l-e~ a _ f . sin(/*) 

— ^7 — ; r- = — / e — ; ^ds, for any r G (0,T , 

Vr(-lnr) vr 2 S v^(-lnr) 1 J ' 

where /* = f*((^s - a 2 )5), i = 1, 2. It is easy to verify that 

lim/* = 0, lim/ 2 * = 0, 

sin (ft) ft a 

hm = 2 ; = hm 



r^o+ ^(-lnr) r-)-o y^(-lnr) 2V/2S ' 

for any s > 0. Using the Lebesgue dominated convergence theorem we finally 
obtain 

E - g Zhu {r) Ea f 00 1 - e~ s , 1 



lim — —. — — = / 77- ds = ■ Ea , 

r^o+ y/r(-lnr) vr 7 (2 s )f ^ 

as claimed. As a consequence of the previous result we can conclude the 
following asymptotic approximation of the formula by Zhu: 

g Zhu (r) esE\l- -j=yft{- In r) J for < r < 1, (5) 

i.e. g Zhu (r) = E (l - ^Vr{- In r)) + o(v/r(- In r)) as r -> 0+ In Fig. ED 

we present a comparison of the analytic solution g Zhu {r) and its asymptotic 
approximation © for r G [0, T] and £ = 100, a = 0.3, r = 0.1, T — 10~ 4 . 



Convexity of the early exercise boundary obtained from Zhu's for- 
mula 

One of the important features of the early exercise boundary for the American 
put option is the convexity of the function g(r) = Sf(T — r) for r G (0, T]. 
The analytic proof of the convexity of g has been presented just recently by 
Chadam et al. in [5]. We also recall that the early exercise boundary is 
log-concave as a function of log of the underlying asset price (cf. Ekstrom 
and Tysk 031112]). 
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Figure 1: Comparison of the analytic solution g Zhu (solid curve) and its 
asymptotic approximation <^ (dashed curve). 
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Figure 2: A graph of the function ft = f^iC'il) 101 various values of the 
parameter 7 (left). A graph of the function £(7) = max^ >0 ft (0 7) (right). 

A relatively simple proof of the convexity of g = g Zhu follows directly 
from the analytic valuation formula Indeed, for any < r < T, we have 
the following expression for the second derivative of the function g Zhu {r): 

In what follows, we shall prove /2(C) = /2(C) 7) £ IP? n ] provided that 7 > 70 
where 70 > is a constant. 

70 = min(7 > | max/ 2 *(C,7) < tt) (6) 

The numerical value of 70 can be estimated as 70 ~ 0.0167821. 

Corollary 2.1 If ^ = 7 > 7o where 70 « 0.0167821 then ft = / S *(C,7) G 
[0,7r] for any ( > 0. yls a consequence, we have ■j^Q Zhu {f) > 0, i.e. the 
function g Zhu (r) as well the early exercise boundary Sf(t) for the American 
put option are convex functions. 
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Remark 2.2 Notice that the condition ^ = 7 > 70 is fulfilled for typical 
market-based choices of the model parameters r and a. For example, if r = 
0.01 (i.e. r = 1% p. a.) then \ > 70 provided that the condition a 2 < 1.19 
(i.e. cr 2 < 119% p. a.) is satisfied. In Figl^we present graphs of the function 
C l— >* /2(C) 7) f or various values of the parameter 7, including the critical 
value 7 = 70 ~ 0.0167821 for which the function £(7) = max^ >0 f^ (C> 7) 
attains the critical value G(7o) = ir. 



2.2 Approximation formula due to Stamicar, Sevcovic 
and Chadam 

In [31] Stamicar, Sevcovic and Chadam derived a single nonlinear inte- 
gral equation for the early exercise boundary position. Based on this inte- 
gral equation the authors derived improved analytical approximation of free 
boundary near the expiry. Asymptotic behavior and justification of the early 
exercise behavior close to expiry have been recently analyzed by Chadam et 
al. in [HI [9] . They proved that the right asymptotic expansion can be ob- 
tained from the nonlinear integral equation developed by Stamicar, Sevcovic 
and Chadam in [3T|. We briefly recall key steps of derivation of the nonlinear 
integral equation for the early exercise boundary position g{r) = Sf(T — r) 
for the free boundary problem ([2]). Let us introduce the following change of 
variables x = In (S/ q(t)) where r = T — t, q(t) = Sf(T — r). Similarly as in 
the case of a call option (see [29]) we define a synthetised portfolio II for the 
put option II(x, r) = V(S, t) — S^(S, t). Then it is easy to verify that II is 
a solution to the following parabolic equation: 

<9IT . N <9II a 2 d 2 Tl „ , _ N 

11(0, r) = E, II(oo, r) = 0, U(x, 0) = 0, x > 0, r e (0, T), (7) 
y^(0,r) = -rE, forre(0,T), 

where a(r) = 44 + r — ^ (see Stamicar ei a/. [HI] , or Sevcovic [221 ED])- 
Applying the Fourier transform one can find the Fourier image of the function 
II in terms of the free boundary position g. The resulting equation for the 
free boundary position reads as ^rfj^(0, r) = —rE, from which the weakly 
singular integral equation for the function g can be found by using the inverse 
Fourier transform (see [21] for details). More precisely, the function g(r) 
fulfills the equation: 



q(t) = Ee~ {r - 



^■)r+<T V2T7?( r ) 
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where the auxiliary function r/(r) is a solution to the following nonlinear 
integral equation 



7}{t) 



\ 



-In 



*V2 



ITT 



a 



forre[0,T]. (9) 



Here the function F v depends on rj via the expression 



F v (r) 



/■tt/2 

2 e 
Jo 
1 



v/2 



COS6 1 



[77 (r) — t)(t sin 2 0) sin 



11^ 



for r G [0, T],# G [0, |]. According to [31], the asymptotic analysis of the 
above integral equation for the unknown function 77 (r) enables us to conclude 
the asymptotic approximation formula for t)(t) as r — > 0. The early exercise 
behavior of q{t) for t — > can be then deduced from the second order 
iteration to the system Q and (11 Op when starting from the initial guess 

2 1 

Vo(t) = (r — corresponding to the constant early exercise boundary 

Sf (t) = E. One can iteratively compute F vo , rji and F m , 7/2 . It turned out 
from calculation performed in [31] that the second consecutive iterate 772 is 
the lowest order (in r) approximation of 77. Namely, 



u 


\2r , 1 

— V27rre rr 




cr 



as r — 7- 0' 



(12) 



Interestingly enough, it has been shown just recently by Chen et al. [8] 
that the early exercise boundary function g is convex (see also [H[9]). More- 
over, the approximation formula (1121) derived by Stamicar, Sevcovic and 
Chadam [21] provides the right asymptotic behavior for r — > + . Further- 
more, Chen and Chadam [9] derived sixth-th order expansion of the function 



17 



51 



287 



199 



24£ 3 64£ 4 120£ 5 32£ e 



(13) 



for £ = In 



87rr 2 



- — > — 00 as r — )• + where 



q{t) = Ee~ a V^). 



(14) 
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Early exercise boundary asymptotic close to expiry 



Similarly as in the case of the analytic approximation formula by Zhu, we 
examine the asymptotic behavior of the function q(t) for r — > where q(t) = 
g (t) is given by the equation: 



Q 



V) = Ee-tr-V^^W , where fj( 



u 


~2r , 1 

— V27rre rr 




a 



(15) 



Employing expression ( II 5p it is straightforward to verify that 

lim v ' = Eg. 

t^o+ ,/tv/ — In r 



Again, as a consequence of the above limit we conclude the following asymp- 
totic approximation of the analytic valuation formula due to Stamicar, Sevcovic 
and Chadam: 



Q SSC 



» &E(l- ay/^V- In t) for < r < 1, (16) 



i.e. q ssc (t) = E(l- ay/r^- In r) + o(y/Ty/-]nr) as r ^ 0+. Notice that 
the asymptotic formula ( 1161) differs from the one obtained from Zhu's formula 
(jlij) by a logarithmic factor v — In r. 



2.3 Approximation formulae by Evans, Kuske and Keller 



In [22] Kuske and Keller proposed another analytic approximation of the 
early exercise boundary for times close to expiration. Then, in the paper 
with Evans [15], they improved and extended the formula for the case of 
dividend-paying asset. 

We begin with the approximation formula by Kuske and Keller [22] . Their 
approximation formula for the position of the early exercise boundary close 
to expiry t — > T reads as follows: 



q (t) ~ £/ l-aV2r 



\ 



In 



2r 9ttt 



as r^0+. (17) 



a V 2 

In [15] Evans, Kuske and Keller derived an improved asymptotic formula: 



Q 



EKK, 



t) w E 1 - <jv2t\ / — In 



"2r 



a 



V2 



TXT 



as t-} . 
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Although, asymptotic formulae (|17p . (|18p by Evans, Kuske and Keller and 
Stamicar, Sevcovic and Chadam (1151) differ in higher order terms of r, it 
holds 

r^0+ ^V-lnr r^0+ ^V-lnT r^0+ y/Ty/-\nT 

It means that approximation formulae due to Evans, Kuske and Keller [22|[T5] 
and Stamicar, Sevcovic and Chadam [31] have the same asymptotic behavior 
close to expiry t ~ T, i.e. < r <C 1. 

3 Numerical methods for calculation of the 
early exercise boundary 

The early exercise boundary function g(r) for the entire time interval r £ 
[0,T], can be approximated by using numerical methods as well. In this 
section we present two approaches: 1) a new local iterative algorithm based 
on the integral equation due to Stamicar, Sevcovic and Chadam [21]; 2) the 
well-known PSOR method (c.f. Kwok [23]). 

3.1 A new numerical algorithm based on a solution to 
the integral equation 

The aim of this section is to introduce a new numerical algorithm for compu- 
tation of the early exercise boundary of the American put option. It is based 
on a solution to the system of implicit equations (TTOj) . (JTTJ derived by 
Stamicar et al. in [31] . The idea of the proposed algorithm is to sequentially 
compute values of the auxiliary function r\ = r/(r) in nodal points 7j £ [0, T}. 
In contrast to global iterative algorithms which iteratively compute the en- 
tire solution q(t),t £ [0,T], (see e.g. Sevcovic [22]) we only need to find a 
root of a real valued function at each nodal point n. This is due to the form 
of functions F v , g v (see fTTOj) and (JTTJ)) whose values at r £ (0, T] depend only 
on the value T)(r) and the history path (r?(^), < £ < r}. 

Our new algorithm for computation of the approximation of the early 
exercise boundary g(r) = Sf(T — r) reads as follows: 

1. Construct a division = r < T\ < ... < r m = T of the interval [0,T]. 
To this end we can employ either equidistant partition r, = (i/m)T, or 
we can use 7$ = (i/m) 2 T in order to adjust the discretization mesh to 
desired behavior (|T6|) of p(r) close to expiry r w 0. We take m ^> 1 
sufficiently large such that —^2Tvr[e rT1 < 1. 
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2. Compute the value of rji ~ rj(ri) from the analytic approximation for- 
mula f[T2l . i.e. 



Vi = ~\ /-In 



2r 



a 



V / 27TTie rTl 



3. for i = 2, ...,m, compute the value 77^ ~ ^fa) as follows: 

3-1 Construct the mapping Q^r^O) = ^ [rjt - fjfa sin 2 9) sin 9] , 
where 77(7^ sin 2 0) is a linear interpolation function between the 
points (Tj,r]j) and (r J+ i,^- + i) if Tj < 7* sin 2 9 < Tj + i for some 
1 < j < z. If < Tj sin 2 < n then 77 (rj sin 2 9) is given by the 
analytic approximation formula (1 12]) . 

3-2 Construct the mapping jF m (ri): 

As for the numerical quadrature of the above integral we can em- 
ploy the composed Newton-Cotes method of the fourth order with, 
at least, 1000 subintervals. 

3-3 Find the root r/i of the equation: 



a \ \ 7r 



The above equation can be solved using either bisection method, 
Newton's method or any other numerical iterative method for find- 
ing roots of real valued functions. In order to speed-up conver- 
gence we can use already constructed value 7?j_i as a starting point 
for iterations at the time level Tj. 

4. Go to step 3 and repeat the calculation of 7/j for the next value of i 
until i < m. 

5. From discrete values rji,i = 1,2, ...,m, we compute the approximation 
Qi of the early exercise boundary position q(ti) as follows: 

We set £0 = E. The entire profile q{t) = Sf(T — t),t G [0, T], is then 
computed as a linear interpolation function between discrete values 
{(Ti,Qi),i = 0,...,m}. 
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3.2 Approximate solution using the PSOR method 

In this section, we present a brief overview how the early exercise boundary 
can be found using a finite difference numerical approximation method ap- 
plied to the variational inequality ([2]). The method consits in computation 
the option price V(S, t) using the so-called projected successive over relax- 
ation (PSOR) method introduced by Ockendon and Elliot in [TJ]. Having 
computed a solution V(S, t) to the variational inequality ^ we can calcu- 
late the early exercise boundary position. Indeed, given a time t, the critical 
stock price Sf(t) is equal to the maximal stock price S = Sf(t) for which the 
option price is equal to the payoff, i.e. 

S f (t) = max{S > | V(S,t) = (E - S) + }. 

Following Kwok [23j . the idea of the PSOR method is to transform (T5]) by in- 
troducing new variables x = ln(S/E), t = T—t, u(x, r) = e ax+ ^ T V(Ee x , T— 
r) where constants a, (3 are defined by a = ^ — \,fl = f + if + ^2- We denote 
u\ ~ u(ih,jk) the finite difference approximation of a solution to the trans- 
formed variational inequality for i = — n, — 1, 0, 1, .., n, j = l,...,m. The 
spatial and time discretization steps h, k > are chosen such that h = L/n, 
k = T/m, respectively. Here T represents expiration time and L is a suffi- 
ciently large bound for the interval x € (—L,L). For practical purposes we 
can take L w 1 . In each time step j = 1,2, ...,m, a linear complementarity 
problem for the finite difference approximation vector u 1 G M 2n+1 is solved by 
using the iterative successive over relaxation (SOR) method where iterates 
are projected to the transformed pay-off diagram. This is done by taking the 
maximum of the transformed pay-off and computed iteration of a solution 
obtained by the SOR successive iteration. For details we refer the reader to 
[231 pp. 212-224]. 

4 Numerical comparison of the early exercise 
boundary approximations 

4.1 Comparison of approximations close to expiry 

This section focuses on numerical comparison of analytic approximations due 
to Stamicar, Sevcovic and Chadam ([8]), Kuske and Keller (fl7|) . Evans, Kuske 
and Keller (ITS]) . Zhu 03) and our new local iterative algorithm from section 
3.1 for the early exercise boundary for times < r = T-t C 1 close to 
expiry. For computational purposes we chose the volatility a = 30%, risk-free 
interest rate r = 10% p. a., and the strike price E = 100$. 
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Figure 3: Comparison of analytic approximation formulae for various matu- 
rities T on a yearly basis. 



Table 1 : Comparison of the early exercise boundary obtained by analytic approx- 
imation formulae and the iterative algorithm to the benchmark PSOR method. 



T 


EKK 


Early 
Zhu 


exercise boundary 

r) = S / (T-r) 
SSCh-A SSCh 


PSOR 


w.r. 
EKK 


Relative error 
to the PSOR method 
Zhu SSCh-A SSCh 


0.000 01 


99.69 


99.51 


99.69 


99.690 


99.7 


0.01% 


0.19% 


0.01% 


0.01% 


0.000 05 


99.37 


99.03 


99.37 


99.358 


99.4 


0.03% 


0.37% 


0.03% 


0.04% 


0.000 1 


99.14 


98.72 


99.15 


99.111 


99.2 


0.06% 


0.49% 


0.06% 


0.09% 


0.000 5 


98.28 


97.57 


98.29 


98.270 


98.31 


0.03% 


0.76% 


0.02% 


0.04% 


0.001 


97.7 


96.83 


97.72 


97.660 


97.73 


0.03% 


0.92% 


0.01% 


0.07% 


0.01 


95.62 


94.27 


95.69 


95.502 


95.6 


0.02% 


1.39% 


0.09% 


0.10% 


0.01 


94.33 


92.73 


94.43 


94.070 


94.18 


0.16% 


1.54% 


0.27% 


0.11% 


0.04 


91.12 


88.66 


91.31 


90.205 


90.3 


0.90% 


1.82% 


1.12% 


0.11% 


0.1 


89.29 


85.25 


89.42 


86.762 


86.94 


2.70% 


1.93% 


2.86% 


0.20% 



Legend: EKK - Evans, Kuske and Keller 1151 . SSCh-A - Stamicar, Sevcovic and Chadam 1311 . SSCh - our 
new local iterative algorithm from section 3.1, ZHU - Zhu 1331 . PSOR - Projected SOR method 1231 . 
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In Fig. [3] we present quantitative comparison of analytic approxima- 
tion formulae by Kuske and Keller (KK), Evans, Kuske and Keller (EKK), 
Stamicar, Sevcovic and Chadam (SSCh-A), Zhu's formula (ZHU). As a nu- 
merical benchmark solution we chose the PSOR method with n = 1000 
spatial grid points and m = 1000 time steps (see section 13.21) . It should 
be obvious that the approximation formulae KK, EKK and SSCh-A exhibit 
similar behavior w.r. to PSOR for the time close to expiry (see Fig. [3] a), 
b)). On the other hand, on a larger time horizon KK, EKK as well as SSCh- 
A become nondecreasing and Zhu's formula (ZHU) better approximates the 
PSOR solution (see Fig. [3jc), d)). It is also worthwile to note that Zhu's 
formula undershoots the early exercise boundary for small values of r when 
compared to KK, EKK, SSCh-A and PSOR. This phenomenon can be easily 
justified by calculating the limit 

E - q ssc {t) 
lo+ E - Q Zhu [r) 



f ^ ^=1- (20) 



In Table Q] we calculated the early exercise boundary position for EKK, 
ZHU, SSCh-A, PSOR methods and our new local iterative algorithm de- 
scribed in section 3.1 which is labeled as SSCh. We also calculated the 
relative error A met ' lod (r) defined as 

A"^ (r) = U s pso R{t £ r) -> T E [0, T], 

where S^ SOR is the early exercise boundary computed by the PSOR method. 
For r ~ 1 minute EKK, SSCh-A and PSOR methods have almost identical 
values (close to 99.4$) whereas Zhu's boundary position has been calculated 
as 99.03$. On the other hand, other approximations (EKK and SSCh-A) dif- 
fers significantly from early exercise boundary obtained by the PSOR method 
as we enlarge time to expiration r > 0.02. The relative error in the early ex- 
ercise boundary position calculated by Zhu's formula w.r. the PSOR method 
is less than 2%. The best approximation of the early exercise boundary has 
been achieved by our local iterative algorithm SSCh. 

In summary, SSCh-A, KK and EKK analytic approximation formulae are 
suitable for approximation of the early exercise boundary close to expiration 
whereas, for a longer time horizon, it is recommended to use the analyti- 
cal approximation formula by Zhu. The new local iterative approximation 
derived in section 3.1 can be used for both small as well as large time horizon. 
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Table 2: Comparison of the early exercise boundary on a long time horizon. 



T 


Early exercise boundary 

Sf(T-r) 
PSOR Zhu SSCh 


Rel. error w.r. to 
PSOR method 
Zhu SSCh 





100. 


100. 


100. 


0% 


0% 


0.02 


92.8672 


90.8575 


92.3461 


2.16% 


0.56% 


0.04 


90.7707 


88.6563 


90.2088 


2.33% 


0.62% 


0.06 


89.3300 


87.2160 


88.7771 


2.37% 


0.62% 


0.08 


88.2350 


86.1300 


87.6695 


2.39% 


0.64% 


0.1 


87.3279 


85.2538 


86.7636 


2.38% 


0.65% 


0.2 


84.2962 


82.3766 


83.7476 


2.28% 


0.65% 


0.4 


81.0179 


79.3593 


80.4793 


2.05% 


0.66% 


0.6 


79.0571 


77.5961 


78.5391 


1.85% 


0.66% 


0.8 


77.6986 


76.3752 


77.1895 


1.7% 


0.66% 


1 


76.6695 


75.4580 


76.1632 


1.58% 


0.66% 


1.5 


74.9137 


73.8879 


74.4094 


1.37% 


0.67% 


2 


73.8107 


72.8731 


73.2722 


1.27% 


0.73% 


3 


72.5786 


71.6205 


71.8735 


1.32% 


0.97% 


4 


72.0121 


70.8778 


71.0464 


1.58% 


1.34% 


5 


71.7966 


70.3925 


70.5100 


1.96% 


1.79% 



Legend: SSCh - our new local iterative algorithm from section 3.1, ZHU - Zhu I33| . PSOR - 1231 . 

4.2 The long term horizon 

In the long term horizon, i.e. r = T — t ~ 1 year or even more, we can 
no longer use the analytical approximations SSCh-A, KK, EKK designed for 
< r < 1 any more. These solutions loose monotonicity for r w 0.1 and 
they become even undefined for large values of r because of the sign change 
in the logarithm. This is why only Zhu's analytical approximation formula 
for the early exercise boundary ((2]) can be used in the long term horizon. 
We compared Zhu's approximation with two numerical methods described in 
section |3j The first method is our new numerical method (labeled by SSCh) 
based on the integral equation fl9]) which was described in section 13.11 The 
second method is the classical PSOR method described in section 13.21 

As a time horizon, we chose the large expiration time T = 5 years. Other 
model parameters are the same as in the previous section, i.e. E = 100$, 
cr = 30%, r = 10% p. a. The computational results are shown in Fig. H] (left) 
and Table [2j We can observe that the shape of all solutions is very similar. In 
Fig. H] (right) we plotted the relative error with respect to the PSOR method, 
which was used as a benchmark. Zhu's analytic approximation formula has 
the relative error between 1 and 2.5% and it attained local minimum around 
t pd 2 years. This is due to the fact that Zhu's method is slightly undershoot- 
ing the solution close to expiry, i.e. for r « 0. The solution computed by our 
new SSCh scheme shows nearly constant error term until r ~ 2.5 years, then 
the relative error starts to grow up. For r « 5 years, the numerical solution 
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T T 

Figure 4: Comparison of the early exercise boundary position in the long time 
horizon. The early exercise boundary position (left). The relative error with 
respect to the PSOR method (right). 



SSCh is approaching Zhu's approximation. This is due to loose of precision 
of the PSOR method itself when the exact early exercise boundary could be 
closer to SSCh and Zhu's approximation than to the PSOR solution. 



5 Comparison of options prices 

In this section we address the question concerning the difference between 
the American put option price and the approximative option price computed 
with an approximation of the early exercise boundary. More precisely, let 
V am (S, t) be the solution to the free boundary problem (J2J) with the early 
exercise boundary profile Sf. Let us consider a given function Sf represent- 
ing an approximation of the early exercise boundary profile Sf. We denote 
by V app the unique solution to the parabolic equation: 



f)V a PP f)V a PP n 2 f) 2 V a PP 

^ + rS ~^S~ + -2 S dS 2 ~~ rVaPP = ^ tE ^ T ^ $f(t)<S, 
V app (+oo, t) = 0, V app (Sf(t), t) — E — Sf(t) , (21) 
V app (S, T) = (E-S) + . 

Notice that we do not require the solution V app to satisfy the C 1 smooth 
pasting contact condition ^gf- (Sf \t) , t) = -1. In fact, V app is a solution 
to the put barrier option (cf. Kwok [23]) with a given down-and-out barrier 
t ^ S app (t). For asset prices < S < Sf{t) we set 

V app {S,t) = E-S. 
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Figure 5: A profile S i— > V am (S, t) of the American option price and its com- 
parison to the option price V app computed with respect to the approximative 
early exercise boundary S apv {t) = g app (T — t). The corresponding European 
style of an option is labeled by V eu . 

A comparison of the profile S t— > V am (S, t) of the American option price and 
the approximative option price V app is shown in Fig. We also plot the 
common lower bound for both put option prices represented by the plain 
vanilla European put option labeled by V eu . 

Knowing the functions t h-> Sf(t) and t H- S apP {t) it is not difficult to 
calculate the difference V am (S, t) — V app (S, t) between option prices. Indeed, 
using the standard transformation (see e.g. Kwok [2"5] ) 

V am {S, t) = Ee ax+l3T u am {x, r), V app {S, t) = Ee ax+)3T u app {x, r), 

where 

i 2 2 

1 t t r u 

x = HS/E),r = T-t, a = ---, f, = ------ 

taking into account the fact that V am (S, t) = E — S for < S < Sf(t) and 
V app (S,t) = E-SiorO<S< S a f pp (t) we can conclude that u° 
solution to the following Cauchy problems: 

for x > ln(g(T)/E), 
re - ax -f3 T for x < ln(^(r)/E), 



du am 


a 2 d 2 u am 


dr 


2 dx 2 


du app 


a 2 d 2 u app 


dr 


2 dx 2 



for x > \n(g app (T)/E), 
re - ax -p T f or x < l n (^pp( r )/ j E) ) 

defined for — oo <x<oo,0<t<T, where g(r) = Sf(T — r), g app (r) = 
S apP (T —t) . Notice that the difference v(x, r) = u am (x, r)—u app (x, r) satisfies 
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v(x, 0) = for each x 6 R. Using Green's representation formula for a 
solution to a linear parabolic equation we obtain, after some calculations, 
the explicit expression for the difference of option prices 

V am (S,t) - V app (S,t) (22) 



rE / / G(x - s,t - Oe a(x ~ s)+P{T ~°ds 

JO J\n{Q a PP{i)/E) 



where G{x,t) = e'* 2 /^ 2 ^ /V2na 2 t is the Green function. The above dif- 
ference in option prices is always nonnegative because the American option 
price is greater or equal to the price of a down-and-out barrier option with 
the prescribed barrier S^ pp (t) = g app (T — t) (see Kwok [23]). 

If we evaluate this difference at the American option early exercise bound- 
ary position Sf(t) then we obtain a slightly simplified expression: 

V am (S f (t), t) - V app (S f (t),t) (23) 

= tE r e -^|iV(7(r,0)-iV(7(r,0)l « 
Jo 

where t = T — t and 

t , l*^ + (^72)(r-0 . t . lng + (r-aV2)(r-0 



Notice that the difference V am (S, t)—V eu (S, t) of the American and European 
style of put options is rather small. Similarly, the difference V am (S, t) — 
V app (S, t) is small. Therefore it is reasonable to calculate the mispricing error 
V am (S, t) — V app (S, t) with respect to the benchmark mispricing difference 
V am (S,t) - V eu (S,t) evaluated at S = S f (t). To this end, let us introduce 
the following relative mispricing error function: 

erriT - t) = f^t) - V« PP {S } {tlt) 

K } V am {S f {t),t)-V eu {S f {t),t)' y } 

The denominator of 024|) can be easily calculated by recalling that 

V am (S f (t),t) = E — S f (t) and V eu {S,t) = Ee^-^ N{-d 2 ) - SN^d^, 

where 

, _ lnf + (r + a 2 /2)(T-t) _ In f + (r - a 2 /2)(T - t) 

d\ — ; , CLo — ; 
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0.001 0.002 0.003 0.004 0.005 0.006 

T 



0.001 0.002 0.003 0.004 0.005 0.006 

T 



Figure 6: Comparison of the early exercise boundary g{j) and the approx- 
imative early exercise boundary g app = g Zhu obtained from Zhu's formula. 
The model parameters were chosen as: E = l,r = 0.1,a = 0.3 for the time 
r = T — t G (0, 0.006) close to expiration. 

(see Kwok [23]). 

In our practical experiment, we evaluated the relative misspricing error 
function err(r) for the approximation of the early exercise boundary obtained 
by Zhu, i.e. we set g app = g Zhu . In Fig. [6] (left) we plotted the relative error 
e(r) in the early exercise boundary position 

SAt) - S Zhu (t) 

between the true early exercise position Sf(t) = g{T — t) and Zhu's approx- 
imation S Zhu (t) = g Zhu (T — t). We can see that the maximal relative error 
in the early exercise boundary position is only 0.32% and it is attained six 
hours prior expiration. 

The relative error function err(r) for times r = T — t close to expiry (less 
than two days) is depicted in Fig. [6] (right). We can see that the error rapidly 
increases when the time t approaches expiration T. For one day to expiration 
(r = 4 x 10~ 3 ) the error is 15%. It increases beyond 70% as t — > T . This is 
due to the fact that Zhu's approximation underestimates the free boundary 
potion as r = T - t 0+ (see $20} ). 

6 Conclusions 

We presented qualitative and quantitative comparison of analytical approxi- 
mations and numerical methods for computation the early exercise boundary 
position of the American put option paying zero dividends. We also pro- 
posed a new local iterative numerical scheme for construction of the entire 
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early exercise boundary which is based on a solution to a nonlinear integral 
equation. We derived asymptotic behavior of approximation formulae for 
the time close to expiry. We proved that the asymptotic formulae by Evans, 
Kuske and Keller (22J US], Stamicar, Sevcovic and Chadam [31] have the 
same asymptotic behavior close to expiry. We also showed that the analytic 
approximation formula by Zhu [33] has a different asymptotic behavior. On 
the other hand, for a long time horizon, Zhu's formula yields quantitatively 
the same results as those of our new local iterative numerical scheme and the 
numerical benchmark PSOR method. 
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comments and suggestions. This research was supported by the bilateral 
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